1 /* Copyright 2002-2016 CS Systèmes d'Information
2 * Licensed to CS Systèmes d'Information (CS) under one or more
3 * contributor license agreements. See the NOTICE file distributed with
4 * this work for additional information regarding copyright ownership.
5 * CS licenses this file to You under the Apache License, Version 2.0
6 * (the "License"); you may not use this file except in compliance with
7 * the License. You may obtain a copy of the License at
8 *
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 */
17 package org.orekit.bodies;
18
19 import org.hipparchus.geometry.euclidean.threed.Vector3D;
20 import org.hipparchus.util.FastMath;
21 import org.orekit.time.AbsoluteDate;
22 import org.orekit.time.TimeScalesFactory;
23 import org.orekit.utils.Constants;
24
25 /** Factory class for IAU poles.
26 * <p>The pole models provided here come from the <a
27 * href="http://astropedia.astrogeology.usgs.gov/alfresco/d/d/workspace/SpacesStore/28fd9e81-1964-44d6-a58b-fbbf61e64e15/WGCCRE2009reprint.pdf">
28 * 2009 report</a> and the <a href="http://astropedia.astrogeology.usgs.gov/alfresco/d/d/workspace/SpacesStore/04d348b0-eb2b-46a2-abe9-6effacb37763/WGCCRE-Erratum-2011reprint.pdf">
29 * 2011 erratum</a> of the IAU/IAG Working Group on Cartographic Coordinates
30 * and Rotational Elements of the Planets and Satellites (WGCCRE). Note that these value
31 * differ from earliest reports (before 2005).
32 *</p>
33 * @author Luc Maisonobe
34 * @since 5.1
35 */
36 class IAUPoleFactory {
37
38 /** Private constructor.
39 * <p>This class is a utility class, it should neither have a public
40 * nor a default constructor. This private constructor prevents
41 * the compiler from generating one automatically.</p>
42 */
43 private IAUPoleFactory() {
44 }
45
46 /** Get an IAU pole.
47 * @param body body for which the pole is requested
48 * @return IAU pole for the body, or dummy GCRF aligned pole
49 * for barycenters
50 */
51 public static IAUPole getIAUPole(final JPLEphemeridesLoader.EphemerisType body) {
52 switch (body) {
53 case SUN:
54 return new IAUPole() {
55
56 /** Serializable UID. */
57 private static final long serialVersionUID = 5715331729495237139L;
58
59 /** {@inheritDoc }*/
60 public Vector3D getPole(final AbsoluteDate date) {
61 return new Vector3D(FastMath.toRadians(286.13),
62 FastMath.toRadians(63.87));
63 }
64
65 /** {@inheritDoc }*/
66 public double getPrimeMeridianAngle(final AbsoluteDate date) {
67 return FastMath.toRadians(84.176 + 14.1844000 * d(date));
68 }
69
70 };
71 case MERCURY:
72 return new IAUPole() {
73
74 /** Serializable UID. */
75 private static final long serialVersionUID = -5769710119654037007L;
76
77 /** {@inheritDoc }*/
78 public Vector3D getPole(final AbsoluteDate date) {
79 final double t = t(date);
80 return new Vector3D(FastMath.toRadians(281.0097 - 0.0328 * t),
81 FastMath.toRadians( 61.4143 - 0.0049 * t));
82 }
83
84 /** {@inheritDoc }*/
85 public double getPrimeMeridianAngle(final AbsoluteDate date) {
86 final double[] m = computeMi(date);
87 return FastMath.toRadians(329.5469 + 6.1385025 * d(date) +
88 0.00993822 * FastMath.sin(m[0]) -
89 0.00104581 * FastMath.sin(m[1]) -
90 0.00010280 * FastMath.sin(m[2]) -
91 0.00002364 * FastMath.sin(m[3]) -
92 0.00000532 * FastMath.sin(m[4]));
93 }
94
95 /** Compute the Mercury angles M<sub>i</sub>.
96 * @param date date
97 * @return array of Mercury angles, with M<sub>i</sub> stored at index i-1
98 */
99 private double[] computeMi(final AbsoluteDate date) {
100 final double d = d(date);
101 return new double[] {
102 FastMath.toRadians(174.791086 + 4.092335 * d), // M1
103 FastMath.toRadians(349.582171 + 8.184670 * d), // M2
104 FastMath.toRadians(164.373257 + 12.277005 * d), // M3
105 FastMath.toRadians(339.164343 + 16.369340 * d), // M4
106 FastMath.toRadians(153.955429 + 20.461675 * d), // M5
107 };
108 }
109 };
110 case VENUS:
111 return new IAUPole() {
112
113 /** Serializable UID. */
114 private static final long serialVersionUID = 7030506277976648896L;
115
116 /** {@inheritDoc }*/
117 public Vector3D getPole(final AbsoluteDate date) {
118 return new Vector3D(FastMath.toRadians(272.76),
119 FastMath.toRadians(67.16));
120 }
121
122 /** {@inheritDoc }*/
123 public double getPrimeMeridianAngle(final AbsoluteDate date) {
124 return FastMath.toRadians(160.20 - 1.4813688 * d(date));
125 }
126
127 };
128 case EARTH:
129 return new IAUPole() {
130
131 /** Serializable UID. */
132 private static final long serialVersionUID = 6912325697192667056L;
133
134 /** {@inheritDoc }*/
135 public Vector3D getPole(final AbsoluteDate date) {
136 final double t = t(date);
137 return new Vector3D(FastMath.toRadians( 0.00 - 0.641 * t),
138 FastMath.toRadians(90.00 - 0.557 * t));
139 }
140
141 /** {@inheritDoc }*/
142 public double getPrimeMeridianAngle(final AbsoluteDate date) {
143 return FastMath.toRadians(190.147 + 360.9856235 * d(date));
144 }
145
146 };
147 case MOON:
148 return new IAUPole() {
149
150
151 /** Serializable UID. */
152 private static final long serialVersionUID = -1310155975084976571L;
153
154 /** {@inheritDoc }*/
155 public Vector3D getPole(final AbsoluteDate date) {
156 final double[] e = computeEi(date);
157 final double t = t(date);
158 return new Vector3D(FastMath.toRadians(269.9949 + 0.0031 * t - 3.8787 * FastMath.sin(e[0]) -
159 0.1204 * FastMath.sin(e[1]) + 0.0700 * FastMath.sin(e[2]) -
160 0.0172 * FastMath.sin(e[3]) + 0.0072 * FastMath.sin(e[5]) -
161 0.0052 * FastMath.sin(e[9]) + 0.0043 * FastMath.sin(e[12])),
162 FastMath.toRadians( 66.5392 + 0.0130 * t + 1.5419 * FastMath.cos(e[0]) +
163 0.0239 * FastMath.cos(e[1]) - 0.0278 * FastMath.cos(e[2]) +
164 0.0068 * FastMath.cos(e[3]) - 0.0029 * FastMath.cos(e[5]) +
165 0.0009 * FastMath.cos(e[6]) + 0.0008 * FastMath.cos(e[9]) -
166 0.0009 * FastMath.cos(e[12])));
167 }
168
169 /** {@inheritDoc }*/
170 public double getPrimeMeridianAngle(final AbsoluteDate date) {
171 final double[] e = computeEi(date);
172 final double d = d(date);
173 return FastMath.toRadians(38.3213 + (13.17635815 - 1.4e-12 * d) * d + 3.5610 * FastMath.sin(e[0]) +
174 0.1208 * FastMath.sin(e[1]) - 0.0642 * FastMath.sin(e[2]) +
175 0.0158 * FastMath.sin(e[3]) + 0.0252 * FastMath.sin(e[4]) -
176 0.0066 * FastMath.sin(e[5]) - 0.0047 * FastMath.sin(e[6]) -
177 0.0046 * FastMath.sin(e[7]) + 0.0028 * FastMath.sin(e[8]) +
178 0.0052 * FastMath.sin(e[9]) + 0.0040 * FastMath.sin(e[10]) +
179 0.0019 * FastMath.sin(e[11]) - 0.0044 * FastMath.sin(e[12]));
180 }
181
182 /** Compute the Moon angles E<sub>i</sub>.
183 * @param date date
184 * @return array of Moon angles, with E<sub>i</sub> stored at index i-1
185 */
186 private double[] computeEi(final AbsoluteDate date) {
187 final double d = d(date);
188 return new double[] {
189 FastMath.toRadians(125.045 - 0.0529921 * d), // E1
190 FastMath.toRadians(250.089 - 0.1059842 * d), // E2
191 FastMath.toRadians(260.008 + 13.0120009 * d), // E3
192 FastMath.toRadians(176.625 + 13.3407154 * d), // E4
193 FastMath.toRadians(357.529 + 0.9856003 * d), // E5
194 FastMath.toRadians(311.589 + 26.4057084 * d), // E6
195 FastMath.toRadians(134.963 + 13.0649930 * d), // E7
196 FastMath.toRadians(276.617 + 0.3287146 * d), // E8
197 FastMath.toRadians( 34.226 + 1.7484877 * d), // E9
198 FastMath.toRadians( 15.134 - 0.1589763 * d), // E10
199 FastMath.toRadians(119.743 + 0.0036096 * d), // E11
200 FastMath.toRadians(239.961 + 0.1643573 * d), // E12
201 FastMath.toRadians( 25.053 + 12.9590088 * d) // E13
202 };
203 }
204
205 };
206 case MARS:
207 return new IAUPole() {
208
209 /** Serializable UID. */
210 private static final long serialVersionUID = 1471983418540015411L;
211
212 /** {@inheritDoc }*/
213 public Vector3D getPole(final AbsoluteDate date) {
214 final double t = t(date);
215 return new Vector3D(FastMath.toRadians(317.68143 - 0.1061 * t),
216 FastMath.toRadians( 52.88650 - 0.0609 * t));
217 }
218
219 /** {@inheritDoc }*/
220 public double getPrimeMeridianAngle(final AbsoluteDate date) {
221 return FastMath.toRadians(176.630 + 350.89198226 * d(date));
222 }
223
224 };
225 case JUPITER:
226 return new IAUPole() {
227
228 /** Serializable UID. */
229 private static final long serialVersionUID = 6959753758673537524L;
230
231 /** {@inheritDoc }*/
232 public Vector3D getPole(final AbsoluteDate date) {
233
234 final double t = t(date);
235 final double ja = FastMath.toRadians( 99.360714 + 4850.4046 * t);
236 final double jb = FastMath.toRadians(175.895369 + 1191.9605 * t);
237 final double jc = FastMath.toRadians(300.323162 + 262.5475 * t);
238 final double jd = FastMath.toRadians(114.012305 + 6070.2476 * t);
239 final double je = FastMath.toRadians( 49.511251 + 64.3000 * t);
240
241 return new Vector3D(FastMath.toRadians(268.056595 - 0.006499 * t +
242 0.000117 * FastMath.sin(ja) +
243 0.000938 * FastMath.sin(jb) +
244 0.001432 * FastMath.sin(jc) +
245 0.000030 * FastMath.sin(jd) +
246 0.002150 * FastMath.sin(je)),
247 FastMath.toRadians( 64.495303 + 0.002413 * t +
248 0.000050 * FastMath.cos(ja) +
249 0.000404 * FastMath.cos(jb) +
250 0.000617 * FastMath.cos(jc) -
251 0.000013 * FastMath.cos(jd) +
252 0.000926 * FastMath.cos(je)));
253 }
254
255 /** {@inheritDoc }*/
256 public double getPrimeMeridianAngle(final AbsoluteDate date) {
257 return FastMath.toRadians(284.95 + 870.5360000 * d(date));
258 }
259
260 };
261 case SATURN:
262 return new IAUPole() {
263
264 /** Serializable UID. */
265 private static final long serialVersionUID = -1082211873912149774L;
266
267 /** {@inheritDoc }*/
268 public Vector3D getPole(final AbsoluteDate date) {
269 final double t = t(date);
270 return new Vector3D(FastMath.toRadians(40.589 - 0.036 * t),
271 FastMath.toRadians(83.537 - 0.004 * t));
272 }
273
274 /** {@inheritDoc }*/
275 public double getPrimeMeridianAngle(final AbsoluteDate date) {
276 return FastMath.toRadians(38.90 + 810.7939024 * d(date));
277 }
278
279 };
280 case URANUS:
281 return new IAUPole() {
282
283 /** Serializable UID. */
284 private static final long serialVersionUID = 362792230470085154L;
285
286 /** {@inheritDoc }*/
287 public Vector3D getPole(final AbsoluteDate date) {
288 return new Vector3D(FastMath.toRadians(257.311),
289 FastMath.toRadians(-15.175));
290 }
291
292 /** {@inheritDoc }*/
293 public double getPrimeMeridianAngle(final AbsoluteDate date) {
294 return FastMath.toRadians(203.81 - 501.1600928 * d(date));
295 }
296
297 };
298 case NEPTUNE:
299 return new IAUPole() {
300
301 /** Serializable UID. */
302 private static final long serialVersionUID = 560614555734665287L;
303
304 /** {@inheritDoc }*/
305 public Vector3D getPole(final AbsoluteDate date) {
306 final double n = FastMath.toRadians(357.85 + 52.316 * t(date));
307 return new Vector3D(FastMath.toRadians(299.36 + 0.70 * FastMath.sin(n)),
308 FastMath.toRadians( 43.46 - 0.51 * FastMath.cos(n)));
309 }
310
311 /** {@inheritDoc }*/
312 public double getPrimeMeridianAngle(final AbsoluteDate date) {
313 final double n = FastMath.toRadians(357.85 + 52.316 * t(date));
314 return FastMath.toRadians(253.18 + 536.3128492 * d(date) - 0.48 * FastMath.sin(n));
315 }
316
317 };
318 case PLUTO:
319 return new IAUPole() {
320
321 /** Serializable UID. */
322 private static final long serialVersionUID = -1277113129327018062L;
323
324 /** {@inheritDoc }*/
325 public Vector3D getPole(final AbsoluteDate date) {
326 return new Vector3D(FastMath.toRadians(132.993),
327 FastMath.toRadians(-6.163));
328 }
329
330 /** {@inheritDoc }*/
331 public double getPrimeMeridianAngle(final AbsoluteDate date) {
332 return FastMath.toRadians(302.695 + 56.3625225 * d(date));
333 }
334
335 };
336 default:
337 return new GCRFAligned();
338 }
339 }
340
341 /** Compute the interval in julian centuries from standard epoch.
342 * @param date date
343 * @return interval between date and standard epoch in julian centuries
344 */
345 private static double t(final AbsoluteDate date) {
346 return date.offsetFrom(AbsoluteDate.J2000_EPOCH, TimeScalesFactory.getTDB()) / Constants.JULIAN_CENTURY;
347 }
348
349 /** Compute the interval in julian days from standard epoch.
350 * @param date date
351 * @return interval between date and standard epoch in julian days
352 */
353 private static double d(final AbsoluteDate date) {
354 return date.offsetFrom(AbsoluteDate.J2000_EPOCH, TimeScalesFactory.getTDB()) / Constants.JULIAN_DAY;
355 }
356
357 /** Default IAUPole implementation for barycenters.
358 * <p>
359 * This implementation defines directions such that the inertially oriented and body
360 * oriented frames are identical and aligned with GCRF. It is used for example
361 * to define the ICRF.
362 * </p>
363 */
364 private static class GCRFAligned implements IAUPole {
365
366 /** Serializable UID. */
367 private static final long serialVersionUID = 20130327L;
368
369 /** {@inheritDoc} */
370 public Vector3D getPole(final AbsoluteDate date) {
371 return Vector3D.PLUS_K;
372 }
373
374 /** {@inheritDoc} */
375 public double getPrimeMeridianAngle(final AbsoluteDate date) {
376 return 0;
377 }
378
379 }
380
381 }